Ultrashort pulse propagation and the Anderson localization 



S. Gentilini 1 ' 2 , A. Fratalocchi 3 ' 1 , L. Angelani 4 G. Ruocco 1 ' 2 , C. Conti 2 
1 Dept. of Physics, University Sapienza, Piazzale Aldo Moro 5, 1-00185, Rome, Italy 
2 Research Center Soft INFM-CNR, c/o Dep. Physics, 
University Sapienza, Piazzale Aldo Moro 5, 1-00185, Rome, Italy 
3 Research Center "Enrico Fermi", Via Panisperna 89/ A, 1-00184, Rome, Italy 
^Research Center SMC INFM-CNR, c/o Dep. Physics, 
University Sapienza, Piazzale Aldo Moro 5, 1-00185, Rome, Italy 

We investigate the dynamics of a 10 fs light pulse propagating in a random medium by the direct 
solution of the 3D Maxwell equations. Our approach employs molecular dynamics to generate a 
distribution of spherical scatterers and a parallel finite- difference time-domain code for the vectorial 
wave propagation. We calculate the disorder- averaged energy velocity and the decay time of the 
transmitted pulse Versus the localization length for an increasing refractive index. 



As originally discussed by Anderson and more re- 
cently in several articles 0, S, i, S, 0, 0, H S, EqL 
[111, QjL [H, [HI, including Bose-Einstein condensation 
[Til , fig , elastic networks [17], and optical beams [l8[, 
three-dimensional (3D) wave localization may occur in 
the presence of structural randomness. An interesting 
issue is the role of localized states in the propagation of 
ultrashort laser pulses in random media [l9[, eventually 
including nonlinear effects [2o| . 

In this Letter we report on ab-initio computational re- 
sults of ultra-short pulses in random media with increas- 
ing refractive index. Our approach combines Molecu- 
lar Dynamics (MD) and parallel Finite-Difference Time- 
Domain (FDTD) codes: the former to provide realistic 
3D distributions of scatterers with quenched disorder, 
the latter to exactly solve the vectorial Maxwell equa- 
tions [20]. With such MD-FDTD technique, we study 
the propagation of classical waves in 3D disordered di- 
electrics for a varying scattering strength, as obtained by 
changing the scatterer refractive index n. As n grows, the 
effective optical path increases and localization occurs. In 
what follows, we characterize the degree of localization 
of the electromagnetic (EM) field by the inverse partic- 
ipation ratio £ of the 3D energy profile and relate it to 
the transmission delay (expressed in terms of the energy 
velocity v e ) and to the decay constant of the trailing edge 
of the transmitted pulse (expressed in terms of an effec- 
tive diffusion coefficient D). As the strong-localization 
regime is attained, D and ve decrease; in addition the 
spectrum of the transmitted pulse displays several nar- 
row peaks. By relating ve to £, we observe a signature 
of a transition, above which localized states are present. 

Our sample is a distribution of 1000 spherical scatter- 
ers obtained by MD simulations. Particle dimensions are 
chosen in order to match typical samples used in exper- 
iments, as e.g. in [2l], [22[. We use a 50/50 mixture of 
particle diameters 310 nm and 248 nm interacting with 
a generalized Lennard- Jones potential 0; at a filling 
fraction (j) = 0.6, this results in a largely disordered and 
tightly packed particle distribution in a cube with edge 
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FIG. 1: (Color online). Transmitted pulse T(£) in linear 
scale for increasing refractive index n; the filled region is the 
adopted input pulse. T(t) is scaled to unitary peak value for 
any n. The inset shows the time delay AT used to evaluate 
the energy propagation velocity ve- 



L = 2.9 jam. The refractive index of the scatterers n is 
chosen in the experimentally accessible range between 1.2 
and 2.8. Several MD runs furnish different configurations 
of the disorder. 

For each realization, we solve the Maxwell equations 
by a parallel FDTD algorithm [23] and study the trans- 
mission of a Gaussian TEMoo linearly y— polarized input 
pulse, with waist wq = 1 /im, impinging on the xy face 
of the cube at normal incidence. The input pulse tempo- 
ral profile is also Gaussian, with duration t = 10 fs and 
carrier wavelength Ao = 532 nm. Numerical results have 
been averaged over five MD configurations of the 1000 
colloidal particles. For each set of MD-FDTD simulations 
we collect: i) the total transmission T(t) by integrating 
the z— component of the Poynting vector with respect to 
the transverse (x, y) output plane, ii) the component 
of the output electric field E y , iii) the corresponding spec- 
trum and iv) the EM energy density £ . In addition, we 
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FIG. 2: (Color online), (a-c) spatial distribution of the energy 
at the output face of the sample (at t = 0.5 ps for a CW 
excitation at A = 532 nm). (d-f) spectrum of the electric field 
E y (thin line) calculated for n — 1.4, 2.2 and n — 2.8; the 
thick line is the spectrum of the incident pulsed beam. 



calculate the distribution of decay times g(r): 



T(t) 



= / e- tlT g(r)dT, 

JO 



(1) 



by a best-fit with a superposition of exponentially de- 
caying functions. The mean value f for the decay-time 
is then used to calculate the "effective" light diffusion 
constant as D = L 2 /7r 2 f (in the broadband short-sample 
regime here considered the diffusion approximation is not 
expected to be strictly valid). The energy propagation ve- 
locity ve is calculated by determining the time AT spent 
by the pulse peak (from the Poynting vector) to travel 
from the input to the output face of the sample (see in- 
set in Fig. [T|), letting ve = L/AT, and averaging over 
disorder realization. 

To characterize the localization length, we calculate the 
inverse participation ratio £: 



J v £ i (x,y,z)dV 



(2) 



being V — L 3 the sample volume. £ is such that if the 
energy profile decays as an exponential with decay con- 
stant it is £ = 2/^; hence it directly measures the 
spatial extension of £(x, y,z). £ is calculated by using a 
continuous wave (CW) beam at A = 532 nm, to avoid 
the simultaneous excitation of several modes. 

Figure [T] shows the input (filled region) and the trans- 
mitted pulses T(t) (solid lines), for increasing values of n; 
the T(t) tail gets longer due to the reduced light diffusiv- 
ity D, while the transmitted pulse slows down (see figure 
[3] below). As shown in Fig. E^i-c, spatial distribution of 
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FIG. 3: (Color online). Phase v p (circles) and energy ve 
(triangles) velocity versus n. 



the energy 6 radically changes from extended (n = 1.4) 
to localized states (n = 2.2 and n = 2.8). Correspond- 
ingly, the spectrum of E y splits into multiple modes (Fig. 
[2ji-f), which implies longer lifetimes for the involved EM 
resonances and a dynamic slowing down (Fig. [T]). 

In Fig. [3] we compare the trend of ^ (triangles) 
versus n with that of the corresponding phase veloc- 
ity v p (circles), which is calculated as v p c± c/n [being 
n = (j) n + (1 — 0), the mean refractive index of the col- 
loidal spheres dispersed in air] . The increase of the degree 
of localization is accompanied by the concurrent swelling 
of the discrepancy between ve and v p: which becomes 
appreciable for n greater than a critical value n c = 1.8. 
Previous experimental investigations have reported sig- 
nificant deviations from an exponential transmission for 
an average index n = 1.55 [22]; in our case we have as the 
critical value for the localization n = 1.48. However fur- 
ther work is required for a strict quantitative comparison 
with experiments. 

Figure [H shows the decay constant distribution g(r) 
of trailing edge of the transmission T(i), as calculated 
from Eq. (pQ). The inset of Fig. [H shows the trailing 
edge of the scaled and temporally shifted transmitted 
pulse in logarithmic scale; the reported linear trends ap- 
parently indicate the absence of any sensible deviation 
from a single exponential, which however becomes evi- 
dent from the spreading of g(r). Note that at the local- 
ization, the width of g(r) is comparable to its mean value 
f . 

Figure [5^ shows the participation ratio £ for different 
values of the refractive index n. As expected the EM 
resonances become more localized as the optical path in 
each scatterer increases; however since £ is a volume av- 
eraged quantity, it does not display a discontinuous trend 
versus n. In Fig. [5b we draw the ratio between phase 
and energy velocity v p /ve and the dynamic diffusion D 
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FIG. 4: (Color online). Distribution of decay-times of the 
trailing edge of the transmission T for various n. The inset 
shows the log-scale plot of normalized transmitted signals T(t) 
(trailing edge). 
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FIG. 5: (Color online). (a) Localization length £ versus refrac- 
tive index n; the error bars and average values as calculated 
over 5 realizations of the disorder. The two insets show the 
snapshot of the energy profile at the output face of the sample 
for two values of n. (b) Parametric plot (versus n) of the ratio 
Vp/ve (left scale) and of D (right scale) versus £. 

versus £. This analysis yields a crossover at £(n c ), where 
the discrepancy between ve and v p is evident. Beyond 
the threshold the reduction of the localization length is 
accompanied by the slowing down of the pulse and a sim- 
ple one-to-one relation between £ and ve is evident. It is 
important to stress that D is expected to vanish at the 
Anderson transition for an infinitely extended structure; 
here we find that, for finite size systems, D (and the en- 
ergy velocity ve) directly measures the degree of spatial 
localization. 

In conclusion we reported on what we believe to be the 
first time-resolved analysis of ultra-short light pulses in 
3D disordered samples. Our approach combines molec- 
ular dynamics and finite- difference time-domain codes, 
thus providing a realistic distribution of the scatterers 
and an exact theory of wave propagation. The distri- 
bution of the decay-time is shown to largely spread at 
the Anderson transition. The plot of the ratio between 
the energy and the phase velocity versus the localization 
length displays a critical character; the more localized are 
the excited EM resonances, the slower is the input pulse 



propagation. The diffusion constant is not vanishing at 
the localization transition and it is a direct measure of 
the spatial extension of the EM field. These findings are 
expected to stimulate novel theoretical works and exper- 
iments in the large community dealing with energy prop- 
agation in the presence of disorder, ranging from optics 
to quantum systems. 
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